function [QW_Mat,Coef_Mat]=Moment2QWCoef(QQ_Unit,CentMomOrderTable,Mean_Mat,CentMoment_Mat,Coef_Mat_Initial)

MultiNum    =   size(CentMoment_Mat,2);
QQNum       =   size(QQ_Unit,1);
QQDim       =   size(QQ_Unit,2);
MomNum      =   size(CentMomOrderTable,1)+QQDim;
QW_Mat      =   zeros(QQNum,MultiNum);
% Moment_Mat  =   reshape(Moment_Vec,[MomNum,MultiNum]);
if nargin<=4
    Coef_Mat_Initial    =   zeros(1+MomNum,MultiNum);
end
Coef_Mat    =   zeros(1+MomNum,MultiNum);
for i=1:MultiNum
    [QW_Mat(:,i),Coef_Mat(:,i)] ...
                =   Moment2Coef_Unit(QQ_Unit,CentMomOrderTable,...
                                     Mean_Mat(:,i),CentMoment_Mat(:,i),...
                                     Coef_Mat_Initial(:,i));
end

return

function [QW,Coef]=Moment2Coef_Unit(QQ_Unit,CentMomOrderTable,Mean,CentMoment,Coef_Initial)

QQ_Dim          =   size(QQ_Unit,2);
QQ_Num          =   size(QQ_Unit,1);
Mean_Target     =   Mean;
CentMoment_Target ...
                =   [zeros(QQ_Dim,1);CentMoment];
CentMoment_Num  =   size(CentMomOrderTable,1);

MaxOrder        =   max(max(CentMomOrderTable));
UnitCentMoment  =   zeros(QQ_Num,MaxOrder+1,QQ_Dim);
CentQQ          =   bsxfun(@minus,QQ_Unit,Mean_Target');
for i=1:QQ_Dim
    UnitCentMoment(:,:,i)   =   bsxfun(@power,CentQQ(:,i),(0:MaxOrder));
end

CentMomentTable =   ones(QQ_Num,CentMoment_Num);
for i=1:CentMoment_Num
    for j=1:QQ_Dim
        Order_Temp      =   CentMomOrderTable(i,j)+1;
        CentMomentTable(:,i) ...
                        =   CentMomentTable(:,i).* ...
                            UnitCentMoment(:,Order_Temp,j);
    end
end
PolyTerm        =   bsxfun(@minus,[CentQQ,CentMomentTable],CentMoment_Target');

% TempFun         =   @(PolyCoef) PolyTerm'*exp(PolyTerm*PolyCoef);
% TempFun_Der     =   @(PolyCoef) PolyTerm'* ...
%                                 bsxfun(@times,exp(PolyTerm*PolyCoef),PolyTerm);
TempQW_Fun      =   @(PolyCoef) exp(PolyTerm*PolyCoef);


if nargin>4
    PolyCoef    =   Coef_Initial(2:end);
else
    PolyCoef    =   zeros(CentMoment_Num+QQ_Dim,1);
end

% PolyCoef    =   fsolve(TempFun,PolyCoef); 

QW          =   TempQW_Fun(PolyCoef);
QW          =   QW/sum(QW);
Temp        =   PolyTerm'*QW;
Diff        =   max(abs(Temp));
ItNum           =   0;
MaxItNum        =   40;
while Diff>1e-6 && ItNum<=MaxItNum
    Temp_Der=   PolyTerm'*bsxfun(@times,QW,PolyTerm);
    PolyCoef=   PolyCoef-Temp_Der\Temp;
    QW      =   TempQW_Fun(PolyCoef);
    QW      =   QW/sum(QW);
    Temp    =   PolyTerm'*QW;
    Diff    =   max(abs(Temp));
    ItNum   =   ItNum+1;
end

M0              =   1/sum(exp(PolyTerm*PolyCoef));

Coef            =   [M0;PolyCoef];
QW              =   M0*exp(PolyTerm*PolyCoef);

return